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Abstract 

Propagating fronts arising from bistable reaction-diffusion equations are a purely deterministic effect. 
Stochastic reaction-diffusion processes also show front propagation which coincides with the deterministic 
effect in the limit of small fluctuations (usually, large populations). However, for larger fluctuations 
propagation can be affected. We give an example, based on the classic spruce-budworm model, where 
the direction of wave propagation, i.e., the relative stability of two phases, can be reversed by fluctuations. 



1 Introduction 

The study of front propagation arising from reaction-diffusion equations is a fundamental problem in 
nonequilibrium physics. In this work we consider cases with two equilibria, a bistable system. If such a de- 
terministic system forms a propagating front, we can say that the "stable" state invades the "metastable" 
state, and the theory [T] gives a simple criterion for which is which. 

However, we can also consider a stochastic system where fluctuations play a role, such that in the 
limit of small noise, we approach the deterministic system. Then the situation is not so simple, and 
fluctuations can play a role in determining the velocity and even the direction of motion of the front, 
because fluctuations can give rise to spontaneous transitions between equilibria. That is, fluctuations can 
reverse the stability of the two equilibria. In this paper we will discuss the mechanism for such reversals, 
and give an explicit example based on the well-known spruce budworm problem [2j El [1] . 

Noise-driven transition between equilibria can be treated using the formalism of rare events. In 
systems with stochastic birth-death processes this approach has recently attracted considerable attention. 
In these systems, fluctuations lead may lead to extinction, the phenomenon that is not described by the 
continuum rate equations [4]. The interplay between stochastic and deterministic effects becomes much 
more intriguing when spatial degrees of freedom are introduced and agents can diffuse on a lattice [5]. 
Here we focus on the role of fluctuations in the phenomenon of front propagation 



2 Deterministic definition of stability 

The standard deterministic treatment of moving fronts in reaction-diffusion systems [T] begins with an 
equation of the form: 

U t = Du xx + f(u). (1) 

Bistability means that 

V{u) = - / dwf(w) (2) 

has two minima, u± , 113 separated by a maximum, 112 . We will refer to V as a potential. 

A moving front occurs when the system is divided into regions where u is at different equilibria. For 
example, u = U3 for x < 0, u — Ui for x > 0. To analyze the dynamics, we seek a traveling front solution 
of Eq. |l]): u(x,t) = it(£ = x — vt), where v is the front velocity. Substituting this into Eq. |l]) we get: 

Du" + vu' + f{u) = 0, (3) 



The prime indicates the derivative with respect to £. Multiplying Eq. Q by u and integrating from 
£ = —oo to oo, we have: 

duf(u)/ / d^D[u'} 2 oc V(ui) ~ V(u 3 ). (4) 

x: 

Note that the velocity will be positive, i.e. the state at 113 will invade wi if Vi — V3 is positive, so that 
«3 corresponds to a lower potential: thus U3 is stable, and ui metastable if V3 < Vi. If V\ — Vz = = 
J™ 1 duf(u) the front will not move (area rule). We refer to this as the "stall point"; note that the stall 
point does not depend on D, but only on V. If Vi — V3 is negative, the front will move towards negative 
x, i.e. ui will invade 113. 

We will need to actually compute v below, i.e. we need the denominator in the Eq. Q. Note that 
Eq. (3} can be thought of as a dynamical system; there are two stable fixed points: u = u\,v! = and 
u — 113, u' = 0. To calculate v we need to find the heteroclinic orbit connecting the two fixed points. We 
do this by using a standard shooting method [6]. 

3 Stochastic definition of stability 

Equations like Eq. [T] often arise as a mean-field description of a stochastic process such as population 
dynamics in spatially extended systems. To be specific, consider a birth-death process for a number, 
n(j), of agents that live on sites Xj with birth rate A(n) and death rate /x(n). Bistability means that 
there are three solutions to A(n) = /J,(n) corresponding to equilibrium population sizes on single sites, 
two stable and one unstable. In the spruce budworm problem, to be treated below, there are two possible 
stable states due to predation: m, the refuge state, and 713, the outbreak state. In our example, below, 
we couple the sites diffusively, i.e. by introducing a rate for transfer to nearest-neighbor sites, e.g., 
(n(j), n(j ± 1)) -> (n(j) - 1, n(j ± 1) + 1). 

It is well-known that the mean-field (large n) limit of the dynamics can be treated by scaling by some 
population scale, A, and thinking of n(j)/A — u(x) as a continuous variable. From the master equation 
we find that the average behavior is given by an equation of the form of Eq. with / = (A — [i)/A. 

However, fluctuations introduce another process, spontaneous transitions between m and 713 on a 
single site. The rate of transitions between the two equilibria involves a different potential 015]: 



<E>(u) = — / dwln 



AW 

* - 2 J (aJ+^w)' (5) 

The last line is the limiting form when //(A + fj,) is small (the Fokker-Planck limit). The transition 
time from state 3 to state 1 is of the form t3_»i = t exp(yl[ ( I > (ti2) — $(113)], where t is a slowly-varying 
prefactor. The transition time depends on the barrier height $2 — $3. A similar expression holds for 
transitions from u\ to U3, ti_»3. 

The smaller time will correspond to the smaller barrier. Based on this, we can give a different 
definition of stability: if i3_»i > ti->3 we say that U3 is stable, and vice-versa. For example, 113 will be 
more stable than ui if the barrier is larger, i.e., if $3 < $1. 

This discussion is for a single site. However, spontaneous transitions can move the front if they occur 
at the interface between the two states. If there are transitions far from the front, and they are not too 
frequent, the "hole" will quickly be filled by diffusion. We have seen examples of this in the simulations 
to be described below. 



4 Competition 

The two potentials, V and $ are similar, but not identical. In particular, there is no reason for them to 
give the same result for stability. Thus a deterministic wave could move from 113 to ui , but noise-induced 
barrier climbing could cause it to move the other way. In this regime, the two effects compete. In the 
next section we give a numerical example where this occurs. Also, fluctuations will shift the stall point. 
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Figure 1: Front profile from simulations of the lattice model, 
width of the interface, which is proportional to y/~D. 



The horizontal box size is taken to be the 



We note that there is a natural length scale for the system given by the width of the wavefront, w. On 
dimensional grounds, w oc yD. This leads us to imagine dividing the system into boxes of width w; see 
Figure [l] Focus on the central box at the front interface. Suppose the continuum front is stalled. If the 
probability for spontaneous jumps 113 — > u\ in the left box is larger than the probability for spontaneous 
jumps iti — > U3 in the right box, the discrete front would move to the left, and the system will not be 
stalled. Thus the observed stall point will be different from that for a purely deterministic front. By the 
same token, if the two transition times are equal, the velocity will be exactly deterministic, but for other 
parameters the velocity will be increased or decreased by spontaneous transitions. 

The relative importance of fluctuations usually depends on the size of the population. For large pop- 
ulations on each site, so that 711,713 3> 1, fluctuations are negligible. As populations decrease stochastic 
effects will start to affect the stall point and v. Finally, for finite systems (perhaps of length < w) and 
very large fluctuations the picture of wave motion will break down altogether, and transitions of the 
system as a whole will be the dominant path for transitions from metastable to unstable. 

We will be interested here in the first two regimes, and the development of waves. All of our qualitative 
notions about the effects of fluctuations on waves will be illustrated by an example in the next section. 

In this work we use periodic boundary conditions. We should note that with absorbing boundary 
conditions, u — at ends of the system, a small system may not be bistable at all [3]. 



5 Waves in the spruce budworm problem 

There is a classic model in the literature, the spruce budworm model EH E], which we will use to 
illustrate the effects we have discussed. The model is based on real experience in forestry where it is 
found that a pest which damages balsam fir trees can exist in two states, the "refuge" state (small 
numbers), and the "outbreak" state, larger numbers. Bistability is thought to be due to the non-linear 
effects of birds on controlling the insect population [2]- 
The spatially extended model [3| is written as follows: 



dn ^d 2 n /_ 71 \ Bn 2 



dt dx 2 V K b J 

The first term on the right hand side of Equation I schematically describes migration, assuming random 
motion. The second term describes the usual logistic growth: insects proliferate with rate rj,, but there 
is some carrying capacity Kb, related to available food, which restricts exponential proliferation. The 
third term represents predation by birds, which saturates at large n: birds are not able to consume more 
than some maximal number of insects per unit time. When n is small, the rate of predation is very 
small, since birds prefer other regions with a larger population of insects. Introducing dimensionless time 
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Figure 2: Velocity, v, as a function of birth parameter r for q = 9, D = 2.5. The solid line is computed 
from Eq. Q and the dotted line from simulations of the discrete model. Positive velocity means that the 
outbreak state takes over the whole system; negative velocity means that the refuge state wins. The dashed 
line is v = 0. The "continuum" stall point r* can be found from the area rule. Note that the "discrete" stall 
point does not coincide with the "continuum" one: r* d > r* . For r between r* and r* d fluctuations reverse 
the velocity of the front. 

t = Bt/A, coordinate x = U B / (AD)x, and dimensionless number of insects u = n/A, we arrive at (bars 
are omitted for clarity) 

du d 2 u 



Here, the proliferation parameter is r — rbA/B and the dimensionless carrying capacity parameter is 
q — Kb I A. In case there is no spatial dependence (a single-site problem), there exists a region in the 
parameter space (r, q), such that there are two stable states: m describes the normal population size and 
us Mi corresponds to the outbreak [HE]. We consider now this region in parameter space and apply 
the formalism outlined above to Eq. Q . Figure [5] shows the front velocity as a function of the birth 
parameter r (solid line). For a specific birth parameter, r* , the two states coexist and the front stalls, 
v = 0. 

Next we analyze the same phenomena as a continuous time Markov process for agents (insects) on 
a lattice. Every site of the one-dimensional lattice can be occupied by any number of insects. At each 
time step, a site, j, is picked at random and then an insect on the site is picked at random. It can either 
jump to a neighboring site (to the right or to the left), proliferate, or die with probabilities related to 
the diffusion, birth and death rates on the site. 

Purth = n/ (r b + M + 2-D) 
Pdeath = fJ,/(r b + n + 2D) 

Pright = Pleft = D / (r b + fl + 2D) 

„ _ r b n{j) Bn(j) 

~i<b~ + A*+nur (8) 

The birth parameter is what we called A above. The death rate per insect, /i, represents the negative 
terms in Eq. Q. After every event, the time is advanced by l/[n t otai(rb + ft + 2D)}, where n to tai is the 
total number of insects. Clearly, Eq. Q is the continuum analog of this discrete lattice model. 

We performed many simulations of the discrete model. To compute the front velocity, we first did 
a time-average for a single run; this was done by plotting the space-time diagram, see Figure [3] and 
computing the slope. Then we averaged the results over many (50 — 100) simulations. Figure [2]presents 
the front velocity as a function of the birth parameter r (dotted line). Strikingly, the discrete stall point 
r* A is not equal to the continuum stall point: r* d > r* . 
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Figure 3: Space-time diagram for the single simulation. Each point on the diagram shows the number of 
particles n in the specific site at some specific time; the color is related to value of n. Initially, a region of 
«3, the outbreak state, was surrounded by the n% state. After the two fronts propagated toward each other, 
the refuge state took over the whole system. The velocity is computed from the slope of the front interface. 
The parameters are q — 9, D — 2.5, r = 0.4605. 

This is precisely the result presented above: since the jump probability from ri\ to 713 is not equal to 
the jump probability from 713 to m, spontaneous jumps contribute to the velocity. Since the probability 
decreases exponentially with the population scale A, the stall point shift r* d — r* should tend to zero 
when the number of particles increases. Simulations of the discrete model show that this is indeed the 
case: r* d tends to its continuum value r* as n increases. 

We suggested above that the "continuum" and "discrete" velocities in Figure [2] should be equal when 
the mean transition time from ni to m, ii_>3 equals the time for the reverse transition, t8_>i. Figure [4] 
shows k — [ti^^/ts^i — 1] versus the birth parameter r. Note that the times are equal at approximately 
f — 0.5075. This suggests that at this value of birth parameter r, there is no stochastic correction to the 
continuum velocity. Our simulations support this prediction: Figure [2] shows that the continuum and 
the discrete velocities are equal at r = f. 

As was mentioned above, the stall point obtained from the continuum theory, r* , does not depend 
on the diffusion coefficient. This is not the case in discrete lattice system: Figure [5] shows the v(r) 
dependence for various diffusion coefficients; r d depends on D. We might expect that the shift r d -— r* 
would tend to zero as the diffusion coefficient increases, and effectively coarse-grains the system. We 
might argue as follows: stochastic transitions contribute to the front velocity when they occur in all the 
sites inside the left or the right box, see Figure [f] The typical front width scales as \[T). The larger the 
box, the smaller the probability for such a collective jump. Therefore, as D increases, the stall point shift 
should tend to zero. However, Figure [5] shows a small discrepancy for large D and no sign of convergence 
to r* = 0.4605. In fact, our simulations (not shown) reveal that the equilibria, ni,ns for the discrete 
model do not approach those for the continuum model for large D, but are uniformly shifted by a small 
amount. We do not understand the large D limit for this system. 

There is another effect that we observed in our simulations. We are dealing with quite small numbers, 
n, so that spontaneous transitions are reasonably common. In this case "islands" of the stable state can 
appear ahead of the front and either disappear, presumably because they are smaller than some critical 
nucleation size, or be enveloped by the advancing front. This has been observed previously in other 
stochastic wave-front problems for similar reasons [9]. 
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Figure 4: The parameter k = (ti^/t^i — l) as a function of r. The circle denotes the value of r (r ~ 0.5075) 
when the two times are equal. The other parameters are q = 9, D = 2.5, A = 7. 
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Figure 5: Velocity of front propagation from the discrete lattice model, for various values of D for q = 9, A 
7. 
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6 Discussion 



It is well known that fluctuations can strongly affect front propagation for situations when a stable state 
invades an unstable state |1UII11||T12] . the "pulled case" . Here we are dealing with "pushed" fronts, and the 
qualitative effects are more subtle, and only occur in a limited region of parameter space. Nevertheless, 
we have shown that they can occur. They should be considered whenever small populations of discrete 
agents are involved in a spatially spreading process. 
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